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Abstract. We have developed a new diagnostic tool for 
the study of gamma-ray emission lines from radioactive 
isotopes (such as 26 A1 and 60 Fe) in conjunction with other 
multi- wavelength observables of Galactic clusters, associ- 
ations, and alike objects. O ur evolutionary synthesis mod- 



els are based on the code of |Ccrviho fc Mas-Hesse (1994) 



which has been updated to include recent stellar evolu- 
tion tracks, new stellar atmospheres for OB and WR stars, 
and nucleosynthetic yields from massive stars during hy- 
drostatic burning phases and explosive SN II and SN lb 
events. 

The temporal evolution of 26 A1 and 60 Fe production, 
the equivalent yield of 26 A1 per ionising 07 V star (Y^ 7 v ), 
and other observables are predicted for a coeval popula- 
tion. The main results are: 

— The emission of the 26 A1 1.809 MeV line is charac- 
terised by four phases: stellar wind dominated phase 
{<, 3 Myr), SN lb dominated phase (~ 3-7 Myr), SN II 
dominated phase (~ 7-37 Myr), and exponential decay 
phase 37 Myr). 

— The equivalent yield V 2 *g 7 v is an extremely sensitive 
age indicator for the stellar population which can be 
used to discriminate between Wolf-Rayet star and 
SN II 26 A1 nucleosynthesis in the association. 

— The ratio of the 60 Fe/ 26 Al emissivity is also an age 
indicator that constrains the contribution of explosive 
nucleosynthesis to the total 26 A1 production. 

We also employed our model to estimate the steady 
state nucleosynthesis of a population of solar metallicity. 
In agreement with other works, we predict the following 
relative contributions to the 26 A1 production: ~ 9% from 
stars before the WR phase, ~ 33% from WR stars, ~ 
14% from SN lb, and - 44% from SN II. For 60 Fe we 
estimate that ~ 39% are produced by SN lb while ~ 61% 
come from SN II. Normalising on the total ionising flux of 
the Galaxy, we predict total production rates of 1.5 Mq 
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Myr -1 and 0.8 Mq Myr -1 for 26 A1 and 60 Fe, respectively. 
This corresponds to 1.5 M of 26 A1 and 1.7 Mq of 60 Fe 
in the present interstellar medium. 

To allow for a fully quantitative analysis of existing 
and future multi-wavelength observations, we propose a 
Bayesian approach that allows the inclusion of IMF rich- 
ness effects and observational uncertainties in the anal- 
ysis. In particular, a Monte Carlo technique is adopted 
to estimate probability distributions for all observables of 
interest. We outline the procedure of exploiting these dis- 
tributions by applying our model to a Active massive star 
association. Applications to existing observations of the 
Cygnus and Vela regions will be discussed in companion 
papers. 

Key words: Stars: abundances — Stars: early-type — Su- 
pernovae:general — Open clusters and associations: gen- 
eral — Gamma-rays: observations 



1. Introduction 

OB associations and young open clusters are the most ac- 
tive nucleosynthetic sites in our Galaxy. The combined ac- 
tivity of stellar winds and core-collapse supernovae ejects 
significant amounts of freshly synthesised nuclei into the 
interstellar medium. Among those, radioactive isotopes, 
such as 26 A1 (r = 1.04 x 10 6 yr) or 60 Fe (t = 2.07 x 10 6 yr), 
may eventually be observed by gamma-ray instruments 
through their characteristic decay-line signatures. Their 
observation presents direct evidence of recent nucleosyn- 
thesis activity, which can be used as a powerful diagnostics 
tool for studies of present galactic activity. 

Galactic 1.809 MeV gamma-ray line emission at- 
tributed to the radioactive decay of 26 A1 has been ob- 
served by numerous gamma-ray telescopes, and the de- 
tailed mapping of the emission distribution by the COMP- 
TEL telescope has clearly identified massive stars as the 
source of this radio-isotope (see Prantzos & Diehl 1996 
for a review). The most convincing evidence for a mas- 
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sive star origin comes from the close resemblance between 
the 1.809 MeV and galactic free-free emission, which links 
26 A1 nucleosynthesis to the O star population (Knodlseder 
et al. 1999a; Knodlseder 1999). Since in general a vari- 
ety of distinct massive star populations of different ages, 
sizes, and metallicities contribute to the observed intensi- 
ties along a line of sight, this indicates that the average 
properties of the populations are related. However, the 
correlation between 1.809 MeV and free-free emission also 
holds for regions far from the Galactic centre where only 
few massive star associations contribute to the observed 
emissions. Examples are the Cygnus and the Vela regions 
where localised 1.809 MeV emission enhancements coin- 
cide spatially with maxima of free-free radiation, show- 
ing the same relative intensities as the Galaxy as a whole 
(Knodlseder et al. 1999a). 

To fully exploit in a quantitative manner such exist- 
ing data and in preparation of the upcoming INTEGRAL 
gamma-ray satellite mission the development of new in- 
terpretation tools is necessary. We here present the first re- 
sults from our modified time-dependent multi-wavelength 
evolutionary synthesis models. A similar model was re- 



cently presented by Pliischke et al. (2000). To model the 



gamma-ray luminosities the nucleosynthetic production 
of the long-lived radio-isotopes 26 A1 and 60 Fe, has been 
included in our multi-wavelength code (Cervino & Mas- 
Hesse 1994, Cervino et al. 2000). These isotopes give rise 
to the 1.809 MeV (for 26 A1) and 1.173 MeV and 1.333 MeV 
(for 60 Fe) gamma-ray lines respectively. The model prop- 
erly accounts for the accumulation of radioactive elements 
and their respective decay times. Results from state-of- 
the-art stellar atmospheres of massive stars are included 
to accurately predict the ionising fluxes from these stars, 
which are at the origin of thermal free-free emission. To- 
gether with the other synthesised observables this provides 
a versatile tool for gamma-ray to radio analysis of massive 
star forming regions. 

For comparisons with individual Galactic star form- 
ing regions (e.g. OB associations, clusters, H 11 regions) 
or ensembles of such objects it is not only imperative to 
model the temporal evolution of their properties. The ef- 
fects of small number statistics of the massive star pop- 
ulation must also be taken into account (Cervino et al. 
2000b). Various studies, including the present one, treat 
such effects by means of Monte Carlo simulations (e.g. 
Cervino & Mas-Hesse 1994, McKee & Williams 1997, Oey 
& Clarke 1998). Finally for a fully quantitative and objec- 
tive confrontation with observables an additional step is 
performed here, to our knowledge for the first time in this 
context. The observational constraints (e.g. a known num- 
ber of stars of given spectral types, derived age, distance 
etc.) and their uncertainties are included in a Bayesian ap- 
proach providing probability distributions for all derived 
properties. 

Section 2 describes the ingredients of our synthesis 
code. The main predictions from our models are presented 



in Section 3. Uncertainties are briefly discussed in Sect. 4. 
Our Bayesian approach to model realistic stellar popula- 
tions is presented in Sect. 5. The main conclusions are 
summarised in Sect. 6. 



2. Evolutionary synthesis model 

2.1. Method 

The starting point for our modelling effort is the evolu- 



tionary synthesis code of Cervino et al. (2000), which pre- 
dicts the time-dependent multi-wavelength energy distri- 
bution of a population of discrete stars from radio wave- 
lengths up to the X-ray domain. For this work, we have 
updated the atmosphere models for massive stars using 
the CoStar models (Schaerer & de Koter 1997), and we 
included at mosphere models for the Wolf-Rayet (WR) 
phase from Schmutz et al. (1992) following the prescrip- 
tions of Schaerer & Vacca (1998). In order to predict 
gamma-ray luminosities, we have included chemical yields 
for the radioactive isotopes 26 A1 and 60 Fe that may either 
be produced during hydrostatic nucleosynthesis in the in- 
teriors of massive stars, or during explosive nucleosynthe- 
sis in supernova explosions (cf. Section 2.3). 

The calculations have been done for two different star 
formation laws in order to explore the extreme cases of 
an instantaneous burst (IB) and of a constant star for- 
mation rate (CSFR). For the IB model, an initial popula- 
tion of coevally formed stars has been created based on a 
Monte Carlo method. Using a power-law initial mass func- 
tion (IMF) of slope r as probability density function]], we 
randomly created initial stellar masses within the interval 
2 — 120 Mq until the total number of stars reaches a pre- 
defined limit. The evolution of each star is then calculated 
using the Geneva evolutionary tracks (see below) in time 
steps of 10 5 years up to an age of typically 50 Myr. At each 
time step the spectral energy distribution and the ejected 
26 A1 and 60 Fe yields are computed for each individual star. 
The evolution of spectral types is also followed in order to 
predict the number of O and WR stars in the population. 
Stars that end their lives during a time step are counted 
as supernova explosions (as far as they are more massive 
than 8 Mq), and are removed from the population for 
the next time step. Summing the contributions from all 
individual star results then in predictions for the entire 
population. 

2.2. Evolutionary tracks 



the population is 
stellar tracks from 
Meynct et al. (1997)| (hereafter MAPP97) for stars w ith 



The evolution 
followed using 



of the stars in 
the non-rotating 



initial mass M ini > 25 M^, |Meynet et al. (1994)| for 



15 < Mini/Ms < 20, and |Schaller et al. (1992)| otherwise 



1 The Salpeter IMF has a slope of T = —1.35 in this pre- 
scription. 
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Solar metallicity tracks are used for this work since we 
are interested in predictions for star clusters and OB as- 
sociations in the solar neighbourhood, but in future we 
plan to extend the calculations also to other metallicities. 
The possible alterations when rotation is taken into ac- 
count in the stellar models are briefly discussed in Sect. 
^. Rotating stellar models will be included when complete 
tracks covering all relevant evolutionary phases will be- 
come available. 

The stellar tracks we used have been calculated for 
enhanced mass loss during the massive star evolution un- 
til the end of the WNL phase. This prescription leads to 
an improved agreement between predictions and several 
observed WR properties; in particular, these models can 
account for the variation of the number ratio of WR to 
O type stars as a function of the metallicity in zones of 
constant star formation rate (Maeder & Meynet 1994). 
The relative populations of WN and WC stars observed in 
young starburst regions are also better reproduced when 
models with high mass loss rates are used (Meynet 1995; 
Schaerer et al 1999). 

The models we are using predict a lower initial mass 
limit of 25 Mq for the formation of WR stars. Uncertain- 
ties due to this mass limit will be discussed in Section [|. In 
order to avoid numerical inconsistencies and an unrealis- 
tic behaviour of interpolated tracks around the WR mass 
limit, we have constructed an artificial track at 25.01 Mq 
going through the WR phases. This track together with 
the published 25 Mq track, which does not enter the WR 
phase, allows smooth interpolations in this mass range. 

2.3. Chemical yields 

The prediction of chemical yields for massive stars de- 
pends critically on the assumed stellar physics, such as 
the treatment of convection, rotation, mass-loss, and the 
final supernova explosion. Being aware of these uncertain- 
ties, we do not intend to predict the yields of 26 A1 and 60 Fe 
to better than within a factor of 2 or so (e.g. Prantzos & 
Dichl 1996; Woosley & Heger 1999), but we merely want 
to identify the main characteristics of nucleosynthesis in a 
massive star population. However, by comparing the mod- 
els to real massive star populations, one can inverse the 
problem and try to learn something about massive star 
nucleosynthesis, and possibly better constrain the theo- 
retical models. In this sense, the employed nucleosynthesis 
yields can be seen as a hypothesis, which can be verified 
by comparison to observations. 

Two different sites of nucleosynthesis must be taken 
into account for the prediction of 26 A1 and 60 Fe yields 
from a massive star association. First, the H-burning in 
the core of the stars may produce appreciable amounts of 
26 A1 that may appear at the stellar surface as an effect of 
both internal mixing and removal of the external layers 
by stellar winds. This 26 A1 can then be ejected into the in- 
terstellar medium by the stellar winds. Second, when the 



star explodes in a supernova event, both 26 A1 produced 
during the post H-burning phases and that synthesised 
at the time of the explosion are then expulsed into the 
interstellar medium. 

To obtain the yields from the population synthesis 
code, one needs a series of different initial mass stel- 
lar models computed with the same physical ingredients. 
Since the Geneva tracks stop before the presupernova 
stage and give no predictions concerning the explosive nu- 
cleosynthesis, we must complement these data with the 
yields of 26 A1 and 60 Fe ejected during the supernova event. 

2.3.1. Stellar winds 



Recent 
sis in 



calculations of hydrostatic 26 A1 nucleosynthe- 
mass losing stars have been undertaken by 
and MAPP97 for non-rotating sin- 



Langer et al. (1995 



gle massive stars. Despite the different assumptions that 
both groups made about mass loss and convection, both 
calculations lead to comparable results (MAPP97). We 
here use the MAPP97 calculations, which have the ad- 
vantage of being consistent with the adopted evolutionary 
tracks which have extensively been compared to observa- 
tions (see Sect. 2.2). The dependence of the total mass 
of 26 A1 ejected by WR stellar winds as a function of the 
initial mass is shown in the left panel of Fig. |l|. We have 
assumed no 26 A1 ejection by stellar winds for the 20 Mq 
stellar model. Log-linear interpolation is used to deter- 
mine the yields at intermediate masses. 

Arnould et al. (1997) suggest that some amount of 60 Fe 
could be ejected by WR stellar winds. Typically for a 60 
Mq star they find ~ 10~ 10 M of 60 Fe ejected during 
the WR phase, which is well below the quantity of 60 Fe 
expelled in the final SN explosion (cf. Fig. ||). The stellar 
wind contribution of 60 Fe can therefore safely be neglected 
for our purposes. 

2.3.2. Type II supernova nucleosynthesis 

At solar metallicity stars in the mass range ~ 8 — 
25 Mq will end their lives by core collapse, giv- 
ing rise to type II supernova (SN II) explosions. 
Explosive 26 A1 and 60 Fe nucleosynthesis during these 
events has been calculated by Woosley & Weaver 
(1995; hereafter WW95), |Thiclcmann ct al. (1996")! and 



[Limongi et al. (2000) . Their models differ in the treat 



ment of the pre-supernova evolution, the prescription of 
convection, the employed nuclear reaction networks, and 
the assumed explosion mechanism. For example, while 
WW95 included hydrostatic 26 A1 production in the pre- 



supernova phase in their models, Thielemann et al. (1996) 
only predict explosive nucleosynthesis yields in their mod- 
els. Additionally, WW95 added neutrino driven spallation 
(the so called ^-process) in their reaction network while 
the others ignore this channel (the ^-process may enhance 
26 A1 production due to additional release of protons). 
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Fig. 1. 26 A1 yields used in this work. Left panel: 26 A1 yields as function of initial stellar mass for Wolf-Rayet stars 
(Meynet et al. 1997, MAPP97, triangles) and type II supernovae (Woosley & Weaver 1995, WW95: diamonds). The 
inset shows 26 A1 yields for Helium stars as function of the initial He mass (Woosley et al. 1995, WLW95). Right panel: 
26 A1 yields as function of initial stellar mass after combining the nucleosynthesis models with the evolutionary tracks 
(see text). For type Ib/c supernovae, 3 link parameters (M^ic- thin solid line, ill/: dashed-dotted, Mqo'- dashed) have 
been explored, which lead to comparable 26 A1 yields. The resulting total yield, obtained using Mco, is shown by the 
thick solid line 




10 20 50 100 10 20 50 100 



M (Mo) M (Mo) 

Fig. 2. 60 Fe yields used in this work. Left panel: 60 Fe yields as function of initial stellar mass for type II supernovae 
(Woosley & Weaver 1995, WW95). The inset shows 60 Fe yields for Helium stars as function of the initial He mass 
(Woosley, Langer & Weaver 1995, WLW95). Right panel: 60 Fe yields as function of initial stellar mass after combining 
the nucleosynthesis models with the evolutionary tracks (see text). Same symbols as in Fig. |[ For M > 25 M Q , the 
final yield corresponds to the yield using Mco as a linking parameter 



Of all these m odels, WW95 predict the highest 26 A1 
yields; for example Thielemann et al. (1996) obtain yields 
that are almost a factor of 10 lower. We here adopt the 
WW95 yields for 26 A1 and 60 Fe, whose dependence on the 
initial stellar mass are shown in the left panels of Figs, [l] 
and H respectively. Note that the WW95 yields do not 
reach down to 8 Mq, the assumed initial mass limit for 
SN II. We assume no 26 A1 production for 7 Mq and per- 
form a log-linear interpolation in the mass interval be- 
tween 8 and 11 Mq. Note that in the case of 60 Fe the 
choice of the inferior mass limit for stars undergoing core- 
collapse has an important impact on the results since the 
stars in this mass range may have a considerable contri- 
bution. 



In order to assign a supernova model (calculated from 
stellar models neglecting mass loss and with different lim- 
its of the convective cores) to the adopted stellar models, 
we use Mco j the mass of the Carbon-Oxygen core at the 
end of C-burning, as linking parameter, as suggested by 
[Maeder (1992) . This procedure is based on the hypothesis 
that the relation between Afco and the explosive nucle- 
osynthetic yields does not much depend on the particular 
set of stellar models. In particular for our case this should 
be a reasonable assumption since at the time of the su- 
pernova explosion the main regions of 26 A1 and 60 Fe pro- 
duction are inside the CO core (cf. WW95). Mco from 
the evolutionary tracks was estimated from the fraction 
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M (Mo) 

Fig. 3. Mho, Mqo, and Mf mass versus initial stellar mass 
as determined from the evolutionary tracks. The mass 
ranges in these three parameters covered by the nucle- 
osynthesis models of WLW95 are shown as vertical lines 
in the left. 

of the convective core before the end of He burning^]. Al- 
though, for the tracks in common with Maeder (1992), 
the derived values are somewhat lower than the ones 
tabulated by Maeder (1992), the yields are only slightly 
modified. For the WW95 models we use the Mqo values 



from Portinari ct al. (1998) calculated by subtracting the 



amount of hydrogen and helium in the WW95 tables from 
the initial mass. The resulting 26 A1 and 60 Fe yields are 
shown in the right panels of Figs. |l| and ^, 

2.3.3. Type Ib/c supernova nucleosynthesis 

For stars that go through the Wolf-Rayet phase, the 
above approach of estimating the explosive nucleosynthe- 
sis yields is not longer valid. The mass-loss will consider- 
ably modify the structure of the star prior to explosion, 
leading eventually to a type lb or type Ic supernova explo- 
sion at the end of its lifetime. After the evaporation of the 
hydrogen envelope, such an object may closely resemble a 
Helium star. 

The only computation of explosive nucleosynthesis 
yields for 26 A1 and 60 Fe in such events comes from Woosley, 
Langer & Weaver (1995; hereafter WLW95) who calcu- 
lated the explosion of mass losing Helium stars of initial He 
masses between 4 — 20 Mq. Again, the models of WLW95 
have to be connected to the evolutionary tracks in our 
evolutionary synthesis model. In principle there are three 
different ways how such a link could be done. First, we 
could use the final mass Mf of the MAPP97 evolutionary 
tracks and link them to the final masses of the WLW95 
models. Second, we could use the mass of the He core A/rc 
at the beginning of core He burning and connect them with 



2 Corresponding to point 42 in the tables of the tracks, where 
the central Helium mass fraction is 0.1. This is justified since 
the subsequent evolution during C-burning should not alter 

M co 



the initial masses of the WLW95 models. And third, as for 
SN of type II, Mco could be used. 

The three possible link parameters Mf, Mjj and Mco 
derived from the evolutionary tracks are shown as function 
of the initial stellar mass in Fig. |^. Since Mh c and Mco 
are not all directly available from the tracks they were es- 
timated from the mass fraction of the convective core close 
to the beginning and end of He-burning respectively^]. The 
Mqo values estimated in this manner are found to be ~ 
20 - 40 % lower than Mco from the stellar structure mod- 
els (Foellmi 1997). As shown by Fig. || not all mass ranges 
overlap with available SN lb models of WLW95. Some link 
parameters are therefore of limited practical use. 

The 26 A1 and 60 Fe yields resulting from the use of the 
different link parameters (excluding extrapolations out- 
side the range covered by the WLW95 models) is shown 
in the right panels of Figs. |l| and || (values at > 25 Mq). 
Given the physical and numerical uncertainties in the link 
between the hydrostatic stellar models and the SN lb cal- 
culations the variations are considered to be small. For 
masses ^ 70 Mq , the link using Mco provides a somewhat 
lower 26 A1 yield than using MH e . However, since these stars 
are relatively rare and since stellar wind ejection exceeds 
the SN lb production by up to an order of magnitude, the 
precise explosive yield in this domain is not crucial. For 
consistency with the type II supernova yields, Mco will 
thus be used for all tracks as the link parameter in the 
remainder of this work. 



Note also that Knodlseder (1999), in his estimation of 
the global galactic 2b Al production rate, used a mass in- 
dependent, constant SN Ib/c yield of 6 x 10~ 5 Mq. This 
prescription is in good agreement with our more refined 
treatment, which predicts yields between (4 — 6) x 10~ 5 

Mq. 

3. Model predictions 

We will now present the temporal evolution of some of 
the key predictions of our model, such as the supernova 
rate, the ionising flux, and the 26 A1 and 60 Fe nucleosyn- 



thesis yields. As indicated earlier (Sect. 2.1) we consider 
two different star formation histories: a coeval population 
(instantaneous burst: IB) and a constant star formation 
rate (CSFR). In the present section an analytic description 
of the IMF is used. Realistic populations of OB associa- 
tions and young open clusters, with a limited number of 
member stars, will be discussed in Section ^| We adopt a 
Salpeter IMF (r = -1.35) over the interval 2 - 120 Mq 
variations of the IMF slope will be discussed in Sect. |3.6| . 
Recall that in both the IB and CSFR cases our normal- 
isation yields absolute quantities given per mass of stars 
formed (IB case), and star formation rate, Mg/yr (CSFR 
case). All other predictions shown here, referring to rela- 
tive quantities, are not affected by the adopted normali- 
sation. 



Points 23 and 42 of the tracks. 
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3.1. Supernova rates 

The predicted supernova rates from our models are shown 
in Fig. ||. For the IB law, the supernova activity starts 
with a sharp peak around ~ 4 Myr which then soon turns 
over into a smoothly declining activity, situated around 1 
SN per Gyr and Mq. The peak is due to the fact that 
stars within the mass interval 60 — 120 Mq all have about 
the same lifetime (3.97, 3.95, and 4.05 Myr for a 120, 85, 
and 60 Mq star, respectively), hence the stars within this 
mass range explode at almost the same moment. The su- 
pernova activity ends around ~ 37 Myr, when all stars 
more massive than 8 Mq have vanished. 



3.2. Ionising flux 

The evolution of the ionising flux Qq, defined as the 
number of photons emitted per second with wavelengths 
shorter than 912 A, is shown in Fig. |5|. After the onset of 
star formation Q decreases with age. In the case of an IB, 
the decline is very rapid, reducing during the hrst 12 Myr 
the ionising flux by a factor of 10 3 . This comes from the 
fact that the bulk of ionising flux is provided by stars more 
massive than ~ 20 Mq , which disappear within only a few 
million years after their formation. In the case of a CSFR 
new massive stars constantly replenish the loss of ionising 
photons from stars which disappear. Since the ionising flux 
increases strongly with mass, Qo reaches equilibrium more 
rapidly than the SN rate. 



. CSFR (yr 1 [Mo/Myr] 1 ) - 



20 30 
Age (Myrs) 



Fig. 4. Predicted temporal evolution of the supernova rate 
for the IB (solid) and CSFR (dashed) star formation laws. 



. IB (photons s Mo ) 

. CSFR (photons s~ 1 [Mo/Myr] - 



Age (Myrs) 



Fig. 5. Predicted temporal evolution of the ionising flux 
for the IB (solid) and CSFR (dashed) star formation laws. 



For the CSFR, the onset of the supernova activity is 
much more smooth, turning quickly into an almost con- 
stant rate of ~ 25 SN per Gyr and M /Myi5 



4 The choice of MQ/Myr instead of the commonly used 
Mq/yt is for illustrative purposes only. 



3.3. 26 Al and 60 Fe ejection rates 

The 26 A1 and 60 Fe ejection rates, j/26 and yeo respectively, 
defined as the mass of radio-isotopes ejected per Myr, 
and normalised to the total mass converted into stars, are 
shown in Fig. |^ for IB models. Several ejection peaks are 
seen in the evolution of the 26 A1 rates. The first peak be- 
tween 2 — 3 Myr, which presents the maximum ejection 
rate during the evolution, comes from the onset of strong 
winds for the most massive stars in the population. The 
next peak at ~ 4 Myr is due to the almost simultaneous 
explosion of all stars in the mass i nter val 60 — 120 M as 
type Ib/c supernovae (see Section 3T and Fig. ^). Around 
~ 5 Myr type Ib/c supernovae start to dominate 26 A1 pro- 
duction, until at 7 Myr the first type II supernovae begin 
to explode. The enhanced 26 A1 production of type II SN 
with respect to type Ib/c leads then to the next peak be- 
tween 7 — 8 Myr. From this age on, the 26 A1 ejection rate is 
dominated by type II events. The time-dependence of the 
ejection rate then reflects the mass dependence of type II 
supernova yields (cf. Fig. and the slow decline of the 
supernova rate (cf. Fig. ||). 

For the 60 Fe ejection rate, a complex structure is found, 
where the first peak again reflects the burst of supernova 
explosions around ~ 4 Myr, and the remaining peaks di- 
rectly reflect the dependency of the 60 Fe yield on initial 
stellar mass (cf. Fig. ||). Hence the temporal structure of 
the 60 Fe rate depends mainly on the details of the explo- 
sive nucleosynthesis models, and following the di scus sion 
about the uncertainties of these models 



2.3 



(Sects. 

the structure should not be regarded as a physical predic- 
tion of our model. In particular, the broad bump between 
13—18 Myr mainly reflects the single high-yield point 
in the WW95 models at M = 13 Mq, and modifications 
in the yields due to changes in the assumptions about the 
stellar physics in the nucleosynthesis models (e.g. Woosley 
& Heger 1999) can easily shift or remove this bump. 

Emissivities (i.e. decay rates) Y(t), in units of M 
per Myr, and normalised to the total mass converted into 
stars, have been obtained by integrating the ejection rates 
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10 20 30 40 50 ^ 10 20 30 40 



Age (Myrs) Age (Myrs) 

Fig. 6. Temporal evolution of the predicted 26 A1 (left panel) and 60 Fe (right panel) ejection rates for instantaneous 
burst models. The immediate total rate is shown by the solid line. Contributions from hydrostatic nucleosynthesis 
(dotted line) and from explosive nucleosynthesis (dashed line) are shown. The dashed-dotted line shows the emissivity 
as defined by Eq. (1) in the text. 



y(t) over the past history including the radioactive decay, 
using 



Y(t) 



y{t')exp(-{t-t')/T)dt\ 



(1) 



where r is the mean life of the radioactive isotope (r = 
1.04 Myr for 26 A1 and r = 2.07 Myr for 60 Fe). As seen in 
Fig. |^ the emissivities are much smoother than the ejection 
rates, which is an obvious result of the convolution opera- 
tion given by Eq. |l|. For 26 A1, the emissivity rises sharply 
between 1 — 3 Myr, which is attributed to stellar wind 
ejecta from massive stars. From ~ 3 Myr on, the emis- 
sivity decreases continuously with a secondary maximum 
around 7 — 8 Myr due to the onset of type II supernovae 
mass ejections. The trend of decreasing 26 A1 rate with in- 
creasing age should be a generic feature for 26 A1 nucleosyn- 
thesis in massive star associations, at least for a Salpeter 
IMF, independently of the uncertainties in the nucleosyn- 
thesis models. The reason is that hydrostatic 26 A1 produc- 
tion, which dominates 26 A1 nucleosynthesis in WR stars 
and SN II, decreases generally with decreasing initial mass 
since the number of seed nuclei becomes smaller. Also, 
the SN rate, which defines the number of ejection events 
within a time interval, decreases with time (cf. Fig. |J). 
However, the level of the first maximum (i.e. the WR- 
peak) with respect to the second maximum (i.e. the SN- 
peak) may depend on details of the nucleosynthesis calcu- 
lations. Also, the age at which the second SN-peak (due 
to SN II) occurs depends on the exact mass limit of WR 
star formation, assumed here to be 25 Mq. E.g. lower val- 
ues of M\vr should shift the SN II peak to older ages. A 
similar temporal behaviour of 26 A1 is found in the models 



of Plfischke et al. (2000) 



For 60 Fe, the emissivity is composed of two broad 
bumps (at ~ 7 and ^17 Myr) that are separated by a local 
minimum around ^13 Myr. After the second bump, the 



emissivity decreases smoothly, followed by the exponen- 
tial decline after ~ 37 Myr. Overall, the 60 Fe production 
stays almost constant between 5 — 18 Myr, followed by a 
slow decline, and due to the uncertainties in the nucleosyn- 
thesis calculations, we should only retain this behaviour 
as characteristic. Given our detailed treatment of yields 
from SN II and SN lb, a larger 60 Fe yield is obtained 



when compared to the models of Pluschke et al. (2000) 
This also affects the predicted 60 Fe/ 2b Al ratio. 

3.4. ea Fe/ 26 Al emissivity ratio 

Part of the 26 A1 and 60 Fe are co-produced in the same 
regions within type II supernovae (Timmes et al. 1995), 
hence the observation of gamma-ray lines from both iso- 
topes can be used as powerful diagnostics tool of nucle- 
osynthesis conditions in such events. For this reason, we in- 
vestigate also the time-dependency of the ratio R between 
60 Fe and 26 A1 emissivities defined as R(t) = Y eo (t) /Y 2G (t) . 
The predicted dependency is shown in Fig. [j] for the case 
of an instantaneous burst, and for a continuous star for- 
mation rate. 

For the IB model, R rises steadily as function of time. 
This is due to the fact that 26 A1 ejection decreases with 
increasing age, while the 60 Fe yields stay roughly constant 
during the active phase of the association. After the last 
supernova exploded, at an age around ~ 37 Myr, the ac- 
cumulated yields decay exponentially, leading to an expo- 
nential rise in R with a time scale of 



1 

1~26 



1 

Ten 



2.09 Myr 



(2) 



(T26 = 1-04 Myr and tqo = 2.07 Myr are the mean lifetimes 
of 26 A1 and 60 Fe, respectively). Note that for ages younger 
than ~ 5 Myr, copious 26 A1 production may appear in an 
association while no 60 Fe has been synthesised yet. This 
is due to the fact that no appreciable amounts of 60 Fe are 
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Fig. 7. 60 Fe over 26 A1 emissivity ratio R for the IB (solid) 
and the CSFR model (dashed). The steady state ratio 
amounts to 0.6. 

supposed to be ejected by stellar winds, and 60 Fe appears 
only when the first supernovae begin to explode. Again, 
this should be a generic feature of the time-evolution of 
a coevally formed population that should be independent 
of details in the nucleosynthesis calculations. 

For a constant star formation rate, R rises more 
smoothly and soon turns into its steady state value around 
0.6. Our result is slig htly in excess of the calculations of 
Timmes et al. (1995) who inferred a value of R = 0.38 ± 
0.27 from a chemical evolution calculation for the Galaxy, 
assuming only SN II nucleosynthesis without mass-loss in 
the initial mass range 11—40 Mq and taking the metallic- 
ity gradient of the Galaxy into account. Adopting similar 
assumptions^ we obtain a ratio of R = 0.43, which is much 
closer to their findings. The main difference between the 
Timmes et al. (1995) and our work lies in the treatment 



of mass loss during the hydrostatic burning phases and 
its effect on the presupernova structure, which leads to a 
reduction of 26 A1 by about 30%, while the 60 Fe nucleosyn- 
thesis remains almost the same. 



3.5. Equivalent 07 V star yields 

COMPTEL observations of the 1.809 MeV gamma-ray 
line suggest that the 26 A1 flux is proportional to the num- 
ber of ionising photons (Knodlseder et al. 1999). In order 
to express the prop ortionality factor in convenient units, 
Knodlseder (1999) introduced the "equivalent 07 V star 



26 A1 yield" , Y^ 7 v , as the mass of 26 A1 that is produced by 
a star of spectral type 07 V, assuming that such a star has 
an ionising flux of logQ^ 7 v = 49.05 ph s _1 . This termi- 
nology follows closely the one employed for the analysis of 
starburst galaxies, where the strength of the ionising flux 
is often expressed in terms of equivalent stars of a given 
subtype (e.g. Vacca 1994). We extend here the definition 



5 Stellar models with no mass loss, 26 A1 produced exclusively 
in SNII with a mass range 11-40 Mq, Salpeter IMF from 



also to 60 Fe, where in analogy Y t 
produced by an 07 V star. 

The predicted equivalent 07 V star yields, calculated 
using 

2Q49.05 



Y 



07 V 



(*) 



Qo(t) 



x Y(t) x t, 



(3) 



0.08 to 40 M . 



are shown in Fig. || for the IB and the CSFR model. Ap- 
parently, Y£q 7 v is an extremely sensitive age indicator for 
a coevally formed population (this would be also true for 
Yg^ 7 v , yet the 60 Fe lines have not been detected so far). 

Within 10 Myr, Y^ r v varies by more than 3 orders 
of magnitude. This strong variation is mainly due to the 
rapid evolution of the ionising flux which drops consid- 
erably as soon as the most massive stars in the associa- 
tion vanished (see Section |3.2| ). In comparison, the most 
common age-indicator in Hn regions, the H/3 equivalent 
width, varies only within 3 orders of magnitude within 20 
Myr (see Cervino & Mas-Hesse 1994 for more details). A 
combination of radio observations providing ionising fluxes 
and 1.809 MeV gamma-ray observations should allow to 
obtain age estimates. In particular, although the strong 
time- variation is driven by the rapid drop in ionising flux, 
adding the gamma-ray observations provides a convenient 
normalisation, making the age estimate independent of 
distance, population richness, interstellar extinction, and 
IMF slope (see also Sect. |3~6| ). 

The evolution of Kg 7 v for the IB model can be split 
into four phases: 

1. the stellar wind phase, lasting from the star formation 
burst up to 3 Myr, and which is characterised by a 
steep rise of Y££ v , 

2. the type Ib/c supernova phase, from 3 — 7 Myr, showing 
a flattening in the slope, which comes from a slight 
decline in 26 A1 production together with a rapid decline 
in the ionising flux, 

3. the type II supernova phase, from 7 — 37 Myr, starting 
with a step around 7 Myr due to the most massive 
SN II, followed by a tail of positive slope since the 
ionising flux drops quicker than 26 A1 production, and 

4. the decay phase, after 37 Myr, which is dominated by 
the exponential decay of 26 A1. 

While this general picture should not depend on details 
of the nucleosynthesis and atmosphere models, the exact 
slopes and time intervals may well change for different in- 
put physics. In particular, phase 3 could already start after 
5 Myr if the minimum mass required to form a Wolf-Rayet 
star would be as high as 40 Mq . Figure || also illustrates 
that the equivalent 26 A1 star yield is an excellent discrimi- 
nator between O star nucleosynthesis (i.e. hydrostatic nu- 
cleosynthesis in WR stars and explosive nucleosynthesis in 
type Ib/c SNe) and SN II nucleosynthesis. While phase 1 
and 2 are characterised by low Y" 2 g 7 v values, ranging from 
zero to ~ 5 x 10 -4 Mq, SN II phase is characterised by 
high equivalent yields well above ~ 10 -3 Mq. Due to the 
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Fig. 8. Equivalent 07 V star 26 A1 yield (left) and 60 Fe yield (right) as predicted by the IB (solid) and the CSFR 
(dashed) model. 



order of magnitude difference, it should be relatively easy 
to use Y^ 7 v to discriminate between both contributions. 

Although the evolution of Y^ 7 v for the IB model fol- 
lows roughly that of Y^ 7 v , we cannot easily identify dis- 
tinct phases as for 26 A1. Due to the complex behaviour 
of the 60 Fe yields, a lot of structure is found in the evolu- 
tion of the 60 Fe ejection rates, but no simple characteristic 
trend. Hence, to first order, the evolution of Y" 6 q 7 v mainly 
reflects the fast decline of the ionising flux with increasing 
age. 



Table 1. Steady-state predictions of the equivalent 07 V 
star yields. 





Source 


Y° 7 v (M ) 


Contribution (%) 


y07 V 
1 26 


MS-winds 
WR-winds 
SN Ib/c 
SN II 


0.46 x 10~ 5 
1.60 x 10~ 5 
0.69 x 10" 5 
2.19 x 10~ 5 


9 

33 
14 
44 




total 


4.94 x 10~ 5 


100 


y07 V 
*60 


SN Ib/c 
SN II 


2.18 x 10" 5 
3.45 x 10~ 5 


39 
61 




total 


5.63 x 10" 5 


100 



The models calculated for a constant star formation 
rate allow us to predict steady state equivalent yields, 
which are reached after ~ 20 Myr (cf. Fig. ||). The re- 
sulting steady state values, split into contributions from 
individual source types, are given in Table [l]. In partic- 
ular, stellar wind contributions have been divided into 
yields ejected prior to (MS-winds) and during (WR- 
winds) the Wolf-Rayet phase, respectively. Apparently, 
1/4 of the hydrostatically produced 26 A1 ejected by stellar 
winds comes from before the WR phase, while the rest is 
ejected when the Hydrogen envelope gets entirely lost in 



a W olf-Rayet phase. A s already pointed out by MAPP97 
and Knodlseder (1999) , stellar wind ejection from massive 
stars provide an important (~ 42 %) contribution to the 
global 26 A1 production. Type II supernovae contribute a 
similar amount, while the rest originates from SN Ib/c 
explosions. The exact repartition on the different source 
classes depends, of course, on the nucleosynthesis mod- 
els, but also on the assumed mass limit for WR star for- 
mation, the slope of the IMF, and finally the metallicity 
(Knodlseder 1999). 

Our model predicts a steady-state equivalent 07 V 
star 26 A1 yield of 4.94 x 10~ 5 Mq, which is lower than 
the observed value integrated over the whole Galaxy of 
(1.0 ± 0.3) x 10~ 4 M (Knodlseder 1999). In view of the 
uncertainties involved in the nucleosynthesis calculations, 
the similarity between model and observation is however 
encouraging. In addition, our models were calculated for 
solar metallicity only, whereas the gamma-ray observa- 
tions average over the entire Galaxy, which shows an av- 
erage metallicity of roughly twice the solar value (Prantzos 
& Diehl 1996). Higher metallicities potentially increase the 
26 A1 production by Wolf-Rayet stars, due to an increase in 
mass-loss and the amount of seed nuclei available for 26 A1 
synthesis (e.g. MAPP97). Hence, including metallicity ef- 
fects in our calculations is expected to raise the Y^ 7 v 
estimate, bringing it even closer to the observed value. 

Using the estimated galactic Lyman continuum lumi- 
nosity of Q = 3.5 x 10 53 photons s -1 (Bennett et al. 1994), 
the number of equivalent 07 V stars can be estimated 
to 31 194, and we can predict galactic nucleosynthesis 
yields from our CSFR model. A similar approach has been 



followed by Knodlseder (1999) using a time-independent 
steady-state model for the Galaxy. In Table ^ we compare 
his findings for solar metallicity and Salpeter IMF with 
mass limits 1 - 120 Mq, to our model (Salpeter IMF with 
mass limits 2 - 120 Mq) and fitting the ionising flux to 
the observed value. Overall, the agreement between the 
models is quite satisfactory. Our models predict a total 
Galactic 60 Fe mass of 1.7 Mq, which due to cancellation 
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Table 2. Galactic yield predictions assuming solar metal- 
licity derived in this work, and given by Knodlseder 
(1999). The star formation rate (SFR) is quoted for the 
mass interval 1 — 120 Mqui Knodlseder (1999) work. 





this work 


Knodlseder (1999) 


Y 26 (M Myr" 1 ) 


1.47 


1.5 ±0.3 


Yeo (M© Myr" 1 ) 


0.84 




M 2(i (Mq) 


1.53 


1.6 ±0.3 


Mao (M ) 


1.74 




SN rate (SN century -1 ) 


2.44 




SFR (M yr" 1 ) 


0.96 


1.2 



of various differences, turns out to be very similar to the 



value of Timmes & Woosley (1997) 



3.6. Dependence on IMF slope 

No general consensus exists about the slope of the IMF in 
young massive star associations and related objects (see 
reviews in Gilmore & Howell 1998). For example from an 
analysis of young open clusters and OB association in the 
Milky Way, Massey et al. (1995) derive an average slope 
of T = -1.1 ± 0.1 for stars with masses > 7Mq. For O 
stars within 2.5 kpc from the Sun parmany et al. (1982) 
find r = —1.6. Based on NIR photometry of the massive 



Cyg OB2 association, Knodlseder (2000) found a compa- 
rable slope of r 
revision 



-1.6 ± 0.1. Finally, in his most recent 
Kroupa (2000")| , obtains T = -1.3 ± 0.7 for stars 
with masses > IMq, taking the scatter introduced by 
Poisson noise and the dynamical evolution of star clusters 
into account. 

Throughout this work a Salpeter IMF slope (r = 
— 1.35) has been used for our "standard" models. The de- 
pendence of our results on T are illustrated subsequently. 




Age (Myrs) 



Fig. 9. 26 A1 emissivity for three different IMF slopes. 

Fig. ^ shows the time-dependent 26 A1 emissivity for 
assuming an IMF slope of V = 0.0, —1.35, and —2.0. 



All three curves have been normalised to the mass trans- 
formed into stars in the mass range 2 - 120 M Q Obvi- 
ously, the structure in the time-evolution remains similar, 
but the importance of stellar wind ejecta with respect to 
supernova ejecta depends strongly on T. For T — 0.0 the 
stellar-wind 26 A1 emissivity peak (at ~ 3 Myr) is almost 
one magnitude larger than for the Salpeter law, leading to 
a burst-like lightcurve that is dominated by stellar wind 
products. In contrast, for T = —2.0 the stellar-wind emis- 
sivity is of the same level as the type II supernova emis- 
sivity, leading to an almost 10 Myrs lasting plateau in the 
lightcurve. 

Interestingly the 60 Fe/ 26 Al emissivity ratio R and the 
equivalent 07 V star 26 A1 yield 1^6 V depend very lit- 
tle on the IMF slope, as shown in Fig. [l(^. This is due to 
the fact that both the nucleosynthetic yield and the ionis- 
ing flux show a similar dependence with initial mass. This 
finding indicates that the equivalent 07 V star 26 A1 yield 
should be fairly reliable age indicator for young massive 
star associations. From Fig. |l^ we estimate a typical age 
uncertainty of ±2 Myr due to IMF variations, which is of 
the same order as typical uncertainties obtained for mas- 
sive star associations by isochrone fitting. Also, the IMF 
variations are smaller than the dispersion introduced by 
statistical fluctuations in a finite sample, as we will demon- 
strate for realistic populations in section 5 (cf. Fig. [ll] right 
panel) . 

4. Model uncertainties 

Before we proceed to the description of our technique 
which will be used in future applications to compare our 
model predictions with observations we shall briefly dis- 
cuss potential uncertainties from stellar models which may 
affect our results. Uncertainties arising from the choice of 
nucleosynthesis models have been discussed earlier (Sect. 
|2.3| ) and shall not be repeated here. 

The main uncertainties related to the stellar models 
included in our synthesis are likely the neglect of stellar 
rotation and the possible importance of massive close bi- 
nary stars as sources of 26 A1 production. 

4-.1. Rotation 

Rotation induces numerous dynamical instabilities in stel- 
lar interiors. The related mixing of the chemical species 
can deeply modify the chemical structure of a star and 
its evolution (see the review by Maeder & Meynet 2000), 
and may in particular have important consequences on the 
26 Al production by WR stars (cf. Arnould et al. 1999). At 
present, very few rotating WR models address this ques- 
tion in detail (see Langer et al. 1995 for some preliminary 
results) , so that it is certainly premature to quantitatively 
assess the possible role rotation plays in that respect. How- 
ever, it seems safe to say that rotation increases the quan- 
tity of 26 Al ejected by WR stellar winds (see Arnould et 
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al. 1999). This statement is corroborated by the numer- 



ical simulations of Langer et al. (1995) which show that 
in case of fast rotation, 2t Al production may be enhanced 
by factors of 2 to 3 with respect to non-rotating stars. 
An increase of 26 A1 by a factor of ~ 1.5 over the non- 
rotating case was found recently by Ringger (2000) for a 
60 Mq model with an average rotational velocity during 
the Main Sequence of 170 km/s. 

Rotation also lowers the minimum initial mass of single 
stars which can go through a WR stage. It is thus likely 
to increase the net amount of 26 A1 ejected by WR stars in 
the ISM. However, since in the present work we use stellar 
tracks with enhanced mass loss rates which, to a certain 
extent, mimic some effects of rotation, it may be expected 
that the inclusion of rotation will not drastically alter the 
present results. 

4-2. Binaries 

The case of binaries is more complicated. First of all the 
stellar physics is even more complex than in single stars, 
which adds numerous potential uncertainties. Second the 
frequency of interacting massive close binaries and their 
exact scenario is poorly known. The importance of bina- 
ries in evolutionary synthesis models is therefore still diffi- 
cult to assess (e.g. Mas-Hesse & Cervino 1999, Vanveberen 
1999). 

For the case of binary systems two effects must be 
taken into account: 1) tidal interactions in close binary 
systems, and 2) mass transfer by Roche Lobe Overflow 
and the posterior evolution of the primary and secondary 
stars. Before any mass transfer will occur tidal effects are 
expected to deform the star and therefore induce insta- 
bilities reminiscent of those induced by rotation. To our 
knowledge, the latter effect has never been studied, even 
if it might have important consequences, as for instance 
by homogenizing the stars and thus inhibiting any mass 
transfer. If mass transfer occurs, the removal of part of 
the envelope of the donor star may favour the appear- 
ance of 26 A1 on the surface in a similar fashion as mass 



loss through stellar winds. More important changes may 
occur for the gainer in systems with primaries Mi ^ 
40 Mq, as e.g. shown by the preliminary studies of 



Braun fc Langer (1995"] and [Langer et al. (1998)[ The lat- 
ter suggest for example a scenario in which mass transfer 
onto the secondary leads to a rejuvenation which alters 
its subsequent evolution. For this case they predict an in- 
crease by 2 — 3 orders of magnitudes of the hydrostati- 
cally produced 26 A1 yield due to a reduction of the de- 
lay between production and ejection of 26 A1, which could 
possibly enhance the total 26 A1 production from type II 
supernovae by a factor of about 2 (Langer, priv. commu- 
nication). 

An attempt to include the contribution of binaries to 
the production of 26 A1 by a stellar population was pre- 
sented by Pltischke et al. (2000)[ Given the limited knowl- 
edge on the evolution and nucleosynthesis from binary 
systems and the numerous potential uncertainties affect- 
ing these predictions, their contribution has deliberately 
been neglected in the present models. As for rotation, if 
all possible effects of binarity produce larger amounts of 
26 A1, our model predictions represent a lower limit for the 
26 A1 production. 



5. Realistic populations 

The model predictions presented in Section have been 
calculated for an infinitely rich association of stars in order 
to avoid effects due to small number statistics, in partic- 
ular at the high-mass end of the IMF. In reality, how- 
ever, galactic associations or open clusters have total stel- 
lar masses between a few 100—1000 Mq (Bruch & Sanders 
1983), limiting the number of O stars, i.e. stars with initial 
masses above ~ 20 Mq, to only a few objects. However, 
since these massive stars provide a considerable fraction of 
the nucleosynthesis yields and ionising power of the associ- 
ation, the early evolution of the system will depend rather 
sensitively on the actual distribution of stellar masses. In 
the following we will describe a Bayesian method that al- 
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lows to quantify the uncertainties introduced by a finite 
sample in our predictions. 

5.1. Probability density functions 

Suppose we want to predict the gamma-ray luminosity 
and ionising flux from an association of age 5 ± 1 Myr, 
which today contains 100 stars within an initial mass 
range of 8 - 25 Mq. At 5 Myr, stars above ~ 40 Mq al- 
ready exploded as supernovae, but the population within 
8 — 25 Mq should still be relatively unaffected by evolu- 
tion. A possible initial population may then be estimated 
by randomly selecting stellar masses from a power-law 
IMF of slope r until the number of stars with masses be- 
tween 8 — 25 Mq amounts to 100 (cf. section |2.1| ). This 
population is then statistically identical to the initial pop- 
ulation of the observed association, although it may differ 
in the exact distribution of stellar masses. Repeating the 
sampling provides several possible initial realisations, and 
the statistical variations in the number and distribution 
of stars among these samples reflects then the ignorance 
about the precise initial conditions in the association. 

The calculation of evolutionary synthesis models for 
each of the samples provides then time-dependent model 
predictions, and the variations among the predictions re- 
flect the uncertainty about the precise distribution of the 
initial stellar population. If the number of random sam- 
ples and corresponding evolution synthesis models is suf- 
ficiently high the probability p(x\t) of observing at an age 
t the value x for quantity X can be reasonably well ap- 
proximated by 

. 1 \ a N(x < X < x + Ax) 
p(x\t)Ax « V ~ '-, (4) 

* samples 

where N(x < X < x + Ax) is the number of samples 
for which X was comprised between x and x 4- Ax, and 
-Samples is the total number of samples. p{x\t) is called 
the probability density function (PDF) of X at age t. 

For illustration, time-dependent probability density 
functions for the 26 A1 emissivity and Y^ 7 v are shown 
for our example association in Fig. O. The PDFs have 
been computed from 1000 Monte Carlo samples assuming 
an IMF slope of T = —1.35. During the early evolution 
(< 3 — 5 Myrs), both quantities are subject to considerable 
uncertainties which can be understood as the combined 
effect of comparable lifetime, 26 A1 yield mass dependence, 
and small number statistics at the high-mass end. As men- 
tioned earlier, stars with initial mass between 60 — 120 Mq 
evolve on comparable timescales, providing contempora- 
neous 26 A1 ejection in this mass range. 26 A1 yields, however, 
vary by almost one magnitude within this mass range, 
making the resulting 26 A1 ejection rates crucially depen- 
dent on the particular spectrum of initial masses. Con- 
sequently, the large statistical uncertainties in the mass 
spectrum at the high-mass end translates into a large un- 



certainty in the 26 A1 emissivity and Y^ 1 v during the early 
evolution of the population. 

During the subsequent evolution, the relative uncer- 
tainty in the derived quantities is roughly constant, with a 
minimum around 5 — 10 Myr followed by a slight rise of the 
uncertainty with time. To understand this behaviour, note 
that the relative uncertainty in the number n of stars that 
contribute to 26 A1 production within a time step is given 
by n -1 / 2 , hence the uncertainty increases with decreasing 
71. Indeed, for ages > 7 Myr the 26 A1 production is only 
due to SN explosions, and the slight decrease in the SN 
rate with time (Fig. [|) is at the origin of the uncertainty 
increase. Note that this feature depends on the slope F 
that is chosen for the stellar population: with a steeper 
IMF the decline in the supernova rate would have been 
reduced (or may even turn into an increase), and conse- 
quently the uncertainty increase would become negligible 
(or might even turn into a decrease of the uncertainties). 

5.2. Inclusion of uncertainties 

In the above example we assumed that the number of 
stars within a given mass interval has been determined 
precisely, and that the slope T of the association is known. 
Both assumptions may not be valid in a realistic case. The 
determination of the association richness is certainly sub- 
ject to some uncertainty, due to membership ambiguities, 
stellar confusion, or invisible members hidden by interstel- 
lar obscuration. 

Uncertainties in both quantities can be easily included 
in the analysis by computing the conditional probability 
density function p{x\t, n,T) for a sufficiently fine grid of 
stellar richness n and IMF slope T. The time-dependent 
association PDF is then calculated using 

p(x\t) — J p(n)p(T)p(x\t,n,T)dndT, (5) 

where p(n) and p(T) are prior probability density func- 
tions quantifying the uncertainties in n and r^. Typically, 
the prior PDF could be a Gaussian if the uncertainty is 
symmetric around a mean value, or a bounded constant if 
the uncertainty is specified as lower and upper limits. 

5.3. Inclusion of prior knowledge 

Finally, to predict the characteristics of an association, the 
information about its actual age and its distance should 
be included in the analysis. In the above example, the age 
was estimated to 5 ± 1 Myr, which again can be expressed 
by a prior PDF p(t). Marginalisation leads then to an age 
independent estimate 

p(x) = J p{t)p(x\t)dt. (6) 

6 This technique of removing irrelevant parameters from the 
problem is called marginalisation. 
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Fig. 12. Posterior PDFs for the 1.809 MeV gamma-ray flux from the decay of radioactive 26 Al(left) and the equivalent 
07 V star 26 A1 yield Y^ 7 v (right). The solid lines show the results for an age of 5 ± 1 Myr, the dashed line has been 
calculated for an age uncertainty of 5 — 15 Myr. A distance of 1.0 ± 0.2 kpc has been assumed for the association. 



Equation g can be seen as a smoothing operation applied 
to the time-dependent PDF over a limited age-window, 
specified by the width of the prior PDF. Since our models 
are calculated for instantaneous star formation, Eq. [| can 
also be interpreted as extending the star formation to a 
period which is defined by the width of the prior PDF. 
In this case, p(t) is the star formation law, and the age 
uncertainty is interpreted as an uncertainty in the star 
formation history. 



Similarly, if the distance s of the association is known 
(to some uncertainty - of course), yields and ionising flux 
can be converted to gamma-ray and radio fluxes, using 

= j p(s)p(<£\s)ds (7) 

where p(s) is the prior PDF, and p(&\s) oc p(x)s~ 2 is the 
distance dependent flux PDF. 

Thus, the final outcome of the evolutionary synthe- 
sis model is not a single value, but a probability density 
function, also called the posterior probability density func- 
tion, that reflects all uncertainties related to the investi- 
gated association. Two examples of such posterior PDFs 
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are shown in Fig. [12] for the 1.809 MeV gamma-ray flux 
and the equivalent 07 V star 26 A1 yield Y^ 7 v . The first 
example (solid lines) has been obtained by assuming an 
age uncertainty of 5 ± 1 Myr, expressed by a Gaussian 
prior PDF with mean of 5 Myr and standard deviation 
of 1 Myr. The second example (dashed lines) has been 
derived for an age uncertainty of 5 — 15 Myr, described 
by a constant bounded prior PDF which is non-zero for 
the interval 5-15 Myr. To obtain the 1.809 MeV flux 
from the 26 A1 emissivity a distance of 1.0 ± 0.2 kpc has 
been assumed, implemented as a prior PDF of Gaussian 
shape with mean of 1 kpc and standard deviation of 0.2 
kpc. Note that Y" 2 ° 7 v is not subject to any distance un- 
certainty since it is defined as a flux (or yield) ratio for 
which the actual distance cancels out. 

In many cases, the resulting posterior PDF has ap- 
proximately a Gaussian shape, but the example of Y^ 7 v 
using an age uncertainty of 5 — 15 Myr illustrates that 
the distribution may be much more complex. The pos- 
terior PDFs may then be used to address various ques- 
tions about the investigated association. For example, the 
probability P that the 1.809 MeV flux is comprised in the 
interval [$ in i n , 'I'max] is simply derived by integrating 



P= p($)d$. 

J$ min 



(8) 



Or, inverting the problem, one may derive the flux in- 
terval [$min)$max] that contains 68% of the probability 
distribution by solving 



*min 



and 



p($)d$ = 0.16 



p($)d$ = 0. 



(9) 



(10) 



(this definition is similar to the la errors often quoted in 
classical statistics). 

The posterior PDFs may also be used to make pre- 
dictions about the detectibility of an association by a 
gamma-ray instrument or radio telescope. If the instru- 
ment sensitivity is given as smallest detectible flux limit 
$3 CT , the probability of detecting the association with the 
instrument is given by 



detect 



p($)c?<f>. 



(11) 



6. Conclusions 



We have constructed a new diagnostic tool for the study 
of the radioactive isotopes of 26 A1 and 60 Fe produced in 
massive star forming regions. The main aim of this work 
was to provide a quantitative model for the analysis of 
multi-wavelength observations of OB associations, open 



clusters and alike objects covering the range from gamma- 
rays (e.g. the 1.809 MeV line of 26 A1, 1.173 and 1.333 Mev 
lines of 60 Fe) to radio, and allowing in a fully quantitative 
manner to account for statistical richness effects of massive 
star populations and other observational uncertainties. 

To achieve this goal we have used the evolutionary 
synthesis models of Gerviho & Mas-Hesse (1994), which 



have been updated to include recent Geneva stellar evo- 
lution tracks, new stellar atmospheres for OB and WR 
stars, and nucleosynthetic yields from massive stars dur- 
ing hydrostatic burning phases and explosive SNII and 
SNIb events (see Sect. 2). In particular proper care was 
taken to combine the stellar models including mass loss 
with appropriate presupernova and SN models. 

The temporal evolution of the ejected quantity of 26 A1 
and 60 Fe produced by a coeval population, other observ- 
ables like the total ionising flux and the supernova rate, 
and derived properties is presented (Sect. 3). This yields 
the following main results: 

— The equivalent 07 V star 26 A1 yield (Y$ 7 v ), defined 
as the stellar yield per ionising flux of an 07 V star 
(Knodlseder 1999), shows a particularly strong time 
dependence where four main phases can be distin- 
guished: stellar wind dominated phase (<^ 3 Myr), SN 
lb dominated phase (~ 3-7 Myr), SN II dominated 
phase (~ 7-37 Myr), and exponential decay phase (<; 
37 Myr). The exact age range of each phase is depen- 
dent on the evolutionary tracks and the lower mass 
limit of WR stars. 

— The use of Y^ 7 v is a powerful tool to constrain the 
evolutionary status of star forming regions. This pa- 
rameter is obtained from a combined observation of 
the 7-ray line and the ionising flux (for example in 
form of thermal free-free radio emission) of an associ- 
ation. 

— 60 Fe production starts with a delay of ~ 2 Myr with 
respect to 26 A1 production. The ratio of the 60 Fe/ 26 Al 
emissivities is also an age indicator that constrains the 
contribution of explosive nucleosynthesis to the total 
26 A1 production. 

Calculations for a steady state population (constant 
star formation; Sect. 3.5) at solar metallicity predict the 
following relative contributions to the 26 A1 production: 
~ 9% from stars before the WR phase, ~ 33% from WR 
stars, - 14% from SN lb, and - 44% from SN II. The large 
contribution from stellar wind ejection (~ 42 %) confirms 
earlier studies of MAPP97 and Knodlseder (1999)| using 
similar yields, who predict contributions of 20-70 % and 
~ 40 % respectively. For 60 Fe we estimate that ~ 39% are 
produced by SN lb while ~ 61% come from SN II. Normal- 
ising on the total ionising flux of the Galaxy, we predict to- 
tal production rates of 1.5M Myr -1 and O.8M Myr -1 
for 26 A1 and 60 Fe, respectively. This corresponds to 1.5 
M Q of 26 A1 and 1.7 M Q of 60 Fe in the present interstellar 
medium. 
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As for other chemical evolution models, our calcu- 
lations depend directly on the adopted nucleosynthetic 
yields, which are affected by considerable uncertainties 
(see e.g. Prantzos 1999). The main uncertainties regard- 
ing 26 A1 and 60 Fe have been discussed in Sects. 2 and 4. 
In fact important new insight especially on the physics 
of supernovae is expected from the study of radioactive 
isotopes such as 60 Fe and 44 Ti, which are synthesised in 
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deep layers close to the so-called mass cut separating the 
outer regions from the remnant. Such studies should also 
benefit from the present models. 

Last, but not least, we have presented a Bayesian ap- 
proach to quantify the predicted observables and their un- 
certainty related to richness effects of the IMF in terms 
of probability density functions (Sect. ||). Subsequently 
these functions can be used in combination with prior 
knowledge on observed objects (e.g. age, distance, and 
their uncertainties) to calculate detection probabilities 
and alike quantities. We have already successfully applied 
our models to existing multi wavelength observations of 
the Cygnus and Vela regions. The results will be published 
in companion papers (Knodlseder et al., in preparation; 
Lavraud et al., in preparation). Our tools will be ideal 
to fully exploit the gamma-ray line observations expected 
from the upcoming INTEGRAL satellite. 
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